Soliton in BCS superfluid Fermi gas 
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We analyze a superfluid state of two species gas of fermions 
trapped in a quasi-lD harmonic potential and interacting via 
attractive s-wave collisions. It is shown that the gap equation 
posesses a self consistent solution with an antisymmetric gap 
function. The gap function has a localized soliton or kink in 
the center of the trap. 
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Over the last decade we experienced rapid development 
in trapping and cooling of dilute atomic gases that made 
condensates of bosonic atoms standard objects of investi- 
gation in laboratories of atomic physics [1] . At the same 
time transition to a superfluid state in cold fermionic 
gases seems to be already attainable [2-4]. The remark- 
able level of manipulation and control of dilute bosonic 
gases allowed for a number of spectacular experiments 
with Bose-Einstein condensates (BEC) [1]. Similar ex- 
perimental activity in fermionic superfluidity may be ex- 
pected in the immediate future. 

In the present Letter we investigate superfluid phase 
of a two-species fermionic gas in the weak coupling BCS 
regime trapped in a quasi one-dimensional (ID) harmonic 
trap. The presence of a trapping potential is a necessary 
condition to have superfluidity in an effective one dimen- 
sional system. We show that the gap function can have a 
soliton. Dark solitons in a ID Bose-Einstein condensate 
have been investigated theoretically and experimentally 
for a few years [5,6]. Atomic density drops to zero in 
the dark soliton. In contrast, the quasi ID Fermi super- 
fluid can support a soliton that does not show up in the 
density distribution but only in the gap function. 

Two species of Fermions ip-( x ) an( i V'+( a; ) with an 
attractive mutual s-wave interactions in a quasi-lD har- 
monic trap are described by the Hamiltonian 



H = 



+ oo 



d.r 



Here Hq = — \-§^i + \x 2 — // is the single particle Hamil- 
tonian in the trap units, g > is an effective ID strength 
of attraction. In the usual BCS mean-field approxima- 
tion with the pairing field A(x) — g{ip-\.{x)ij)-(x)) and 
the Hartree-Fock potential W(x) = —g(ip±(x)'tl!±(x)) the 
stationary problem at zero temperature is equivalent to 
solving the Bogoliubov-de Gennes equations [7] 

+H u m + W(x)u m + A(x)v m , 
u m v m = -T-toVm - W(x)v m + A*(x)u m , (1) 



together with the self-consistency conditions 

oo 

oo 

W{x) = -g J2 K(z)| 2 ■ 



(2) 
(3) 



m — 1 



Unlike in 2D or 3D [8,9], in ID both sums are convergent. 
Similar Bogoliubov-de Gennes equations with soliton so- 
lutions appear in the theory of conducting polymers [10] 
where the order parameter A is real and lives on a dicrete 
ID lattice. 

In numerical calculations the infinite sums must be re- 
placed by finite sums up to a certain cut-off A, 



A A (x) = g Y u m {x)v* m (x) 



(4) 



m— 1 



For a sufficiently large cut-off A A (x) converges to A (a;). 
In principle for fixed (u rn ,v rn ) the difference is A (a;) — 
A A (x) = ©(A" 1 / 2 ). Unfortunately, we found that the 
small error 0(A -1 / 2 ) was amplified by the nonlinear set 
of Eqs.(l)-(3) and the convergence was much slower. This 
is why we use a different prescription 



A(x) = A A (x) +5A A (x) 



(5) 



Here SA A (x) w jE^a+i^W"^^) is an approxi- 
mate analytical expression for the part of the infinite sum 
that is normally ignored. Again for A large enough we re- 
cover the infinite sum, but now the convergence is much 
faster than for the simple cut-off in Eq.(4). 

We use the local density approximation (LDA) to eval- 
uate the high energy part of the gap function 



8A A (x) = g 



dk A(x) 



2tt 2E(k,x) 



(6) 



In the LDA the superfluid is considered locally uniform 
so that locally one can use the BCS theory for a uniform 
system. In equation (6) the E(k, x) — \J A 2 (x) + e 2 (k, x) 
is a local quasiparticle energy with e(k, x) = 4- + %- + 
W(x) — fi. A local cut-off momentum k A (x) correspond- 
ing to the cut-off energy + /i is given by k\(x) = 
2u> a + 2/i — x 2 — 2W(x) when positive and zero otherwise. 
A combination of Eqs.(5) and (6) gives to leading order 
in A{x) 
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A 

A(x) w 3A (x) u m (x)v* m (x) , 

m— 1 

i def. i _ j_ r 00 rffc 

Sr A (z) .9 2n J kA(x) 2E(k,x) 



(7) 
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fcA(x) + fc F (x) 

5 2nkp(x) \k\(x) — ftp (a:) 



In 



Here Uf{x) = y/2/j, — x 2 — 2W(x) is a local Fermi mo- 
mentum. The LDA applied to the highest quasiparti- 
cle states results in an ^-dependent attraction strength 
g\(x). The enhanced gh{x) > g approximately compen- 
sates for the high energy quasiparticle states missing in 
the finite sum. We use the LDA to evaluate only the sum 
SAa(x) over the high energy states. Although accuracy 
of LDA as such may be disputable, there is little doubt 
that it is a very accurate approximation in the high en- 
ergy regime where k\(x) is much greater than any rate 
of variation of A(x). 
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FIG. 1. Gap functions A(x) in the antisymmetric state 
(solid lines) and approximate LDA gap functions Alda(;e) 
in the symmetric ground state (dashed lines) for g = 6 [panel 
(a)], g — 8 [panel (b)], g = 10 [panel (c)], and g = 12 [panel 
(d)]. The chemical potential is /i = 50 which gives an average 
total number of atoms ranging from 2 x 70 for g = 6 to 2 x 82 
for g = 12. 

In this paper we consider an antisymmetric A(x) = 
—A (—a;) with a symmetric density of atomic cloud 
W(x) = W(— x) in the weak coupling BCS regime. The 
solid lines in Fig.l show the antisymmetric gap functions 
for four values of the attraction strength g. Panels (b)- 
(d) show very clear characteristic soliton (kink) patterns, 
but in panel (a) the width of the soliton is comparable 
to the size of the atomic cloud and we see just a plain 
antisymmetric state. 

The dashed lines in Fig.l show symmetric ground state 
A(x)'s obtained with the LDA. To get this approximate 



gap function we set A = on the left hand side of 
Eq.(6) which replaces the #Aa(x) by an approximate 
A (a;). Then we solve Eq.(6) to leading order in A(x): 



A LDA (x) = 4 k 2 F {x) e -»M*)/fl 



(8) 



wherever fcp(x) is real and otherwise. The agreement 
between the symmetric Alda and modulus of antisym- 
metric gap functions is reasonable even for the small 
numbers of atoms considered in Fig.l. The solitons in 
panels (b)-(d) can be considered as localized defects in a 
locally uniform superfluid. 

The gap function A(x) has a soliton but, unlike for 
dark solitons in BEC [5,6], corresponding atomic den- 
sity shows hardly any symptoms of the soliton at all, see 
Fig. 2. The soliton in the weakly interacting BCS super- 
fluid is not a dark soliton. 
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FIG. 2. Atomic densities p(x) = ( , 4>±(x) , 4>±(x)) corre- 
sponding to the antisymmetric states in Fig.l. The fast os- 
cillations of p(x) in the solitons average to zero on length 
scales longer than fe^ 1 . The dashed plots show approximate 
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To demonstrate convergence of the scheme based on 
Eq. (7) we show A(x) for g = 12 and increasing value of A 
in Fig. 3. In this figure we also show results corresponding 
to Eq. (4) where the contributions from the high energy 
states are neglected. We see that the results based on the 
scheme Eq. (7) are closer to the true solution, even for 
A = 400, than those based on the simple cut-off scheme 
(4) with A = 1000. 

We need some analytical insights to extrapolate the 
numerical results to a large number o atoms N. The 
most intriguing issue is the width of the soliton. A sim- 
ple estimate of the healing length £ can be obtained in 
the LDA. The local quasiparticle energy E{k,x) can be 
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expanded in small fluctuations 6k = k — kp around Fermi 
momentum 
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(9) 



We see that the quasiparticle spectrum is affected by the 
gap only in the range of 5k w A/kp. This is the disper- 
sion of momenta available to construct a localized soliton 
and so the healing length must scale as 



2fi 



(10) 



with a constant a = 0(1). To estimate the constant 
we fit the solitons in panels (b)-(d) in Fig.l with a trial 
function 7Alda(x) tanh(x/£), where 7 and £ are fitting 
parameters. The best fit to the central part of panel (d) 
is shown in Fig. 4. From £'s we obtain the best fits for a's 
and show them in Table I — it is fair to say that awl. 




FIG. 3. Gap functions calculated in different approxima- 
tions for chemical potential /u = 50 and attraction strength 
g = 12. Dotted line corresponds to Eq. (4) with A = 1000. 
Dashed, dotted-dashed and solid lines are the results that base 
on Eq. (7) with A equal to 400, 700 and 1000, respectively. 




FIG. 4. Solid line indicates the gap function shown in panel 
(d) of Fig. 1. Dashed line is a fit to the solid line with a trial 
function 7Alda(s) tanh(a:/£) perfomed in the range of space 
presented in the figure only. The best fit was obtained for 
a = 1.07, Eq. (10). 

Another interesting question is when the soliton is 
a nice localized kink as in panels (b)-(d) of Fig.l and 
when its width is comparable or exceeds the width of the 
atomic cloud as in panel (a). A simple estimate of the 
"kinkness" of the antisymmetric state is the ratio of the 
healing length to the Thomas-Fermi width of the cloud 
x TF ~ V^, 



XTF 



1 

A 



(11) 



This factor can be made small for large N when a large 
A 3> 1 is easily compatible with the weak coupling con- 
dition A <C fi necessary for the LDA to be accurate. 

Once we have a localized soliton another interesting 
question arises. A localized soliton should behave like a 
non-dispersive point particle. As the atomic cloud in the 
trap is non-uniform this particle feels an effective poten- 
tial V e s(q) with soliton position q. The question is what 
is the shape of this potential. This is an important issue 
related to stability because, for example, the dark soli- 
tons in BEC experience an inverted harmonic potential 
[5]. Their small fluctuation spectrum has an anomalous 
mode with negative frequency which might imply their 
instability. Fortunately, it does not for technical reasons 
specific to the dark solitons in harmonic traps. To get 
an insight into V e e(q) we make a simple estimate of the 
soliton energy. According to Ref. [7] the difference be- 
tween the energy density of the superfluid state with a 
gap A^0 and the normal state with A = is given by 



Se 



-f 

Jo 



dA A 



d g~ 



dA' 



(12) 



When combined with the LDA dependence of the gap 
on the attraction strength g in Eq.(8), this formula gives 



3 



Se = —A 2 /2irkF- Now, in zero order approximation the 
soliton is an area of width £ where this energy gain is 
missing, hence the energy of the soliton is roughly — Se x £ 
or 

V cS (q) ~ |A L DA(z = g)| . (13) 

We come to a remarkably simple conclusion that the ef- 
fective potential is proportional to the dashed plots in 
Fig.l. Small fluctuations of the soliton close to the min- 
imum of the potential in the center of the trap have pos- 
itive frequency. There is no anomalous mode that might 
imply thermodynamic instability. In order to deexcite 
the superfluid from the soliton state to the symmetric 
ground state the soliton has to be pushed beyond the 
potential barrier at the edge of the atomic cloud. In 
this sense solitons in the BCS Fermi superfluid are more 
robust than dark solitons in BEC. Hopefully the dark 
soliton experiments in BEC will find their counterpart in 
BCS Fermi superfluid. 

As was shown in Ref. [11], in a 3D Fermi superfluid 
in an isotropic harmonic trap there is no energy gap in 
the quasiparticle energy spectrum even for a large gap 
function A(x) in the center of the trap. Thus in 3D 
the spectroscopic detection of the BCS state requires the 
probing laser beam to be focused on the locally uniform 
cenral part of the trap. However, the detection schemes 
based on the measurement of the energy gap [12] can 
be very useful in ID. In ID an nonzero A(x) directly 
translates into a gap in the quasiparticle spectrum. 

What is more the same schemes can be used to detect 
the soliton in the superfluid. In Table II we list the lowest 
quasiparticle energies for the solitons in panels (c)-(d) 
of Fig.l. The energies to m >i start roughly at the gap 
close to the center of the trap |A(x s=s 0)|. This is the 
usual gap in the BCS state. However, the lowest u)\ is an 
order of magnitude lower. This is a quasiparticle bound 
state localized inside the soliton where the gap function 
is close to zero. This mode is a soliton counterpart of the 
Caroli-de Gennes bound states inside a vortex core [13]. 
The bound state manifests itself in the density plots in 
Fig. 2 by the fast density oscillations inside the soliton. 
The spectroscopic detection schemes will see the energy 
of the bound state in the soliton and they can be used 
both to detect the BCS state in ID and to distinguish 
the ground state from the soliton state. 

In conclusion, we have explored physics of soliton in su- 
perfluid atomic Fermi gases. Well defined kink-like soli- 
ton in the gap function can be obtained in a quasi- ID 
trap. We have estimated its width and energy and ar- 
gued that small fluctuations of a soliton around the cen- 
ter of the trap have positive frequency. In this respect 
the BCS soliton is more robust than the dark soliton in 
a Bose-Einstein condensate. The soliton in superfluid 
fcrmi gases can be detected by spectroscopic measure- 
ments where the energy of the quasiparticle bound state 
localized in the soliton may be observed. 



TABLE I. Values of a parameter, Eq. (10), obtained from 
fitting of a trial function, 7Alda(s;) tanh(cc/£), to the gap 
functions shown in panels (b)-(d) of Fig. 1. 



Panel in Fig. 1 (b) (c) (d)_ 

a 1.07 1.05 1.07 



TABLE II. The lowest quasiparticle energies for the soli- 
tons in panels (c)-(d) of Fig.l. 

Panel in Fig. 1 (c) (d)_ 

wi 1.9 4.0 

LO'2 11.5 20.1 

lj 3 11-7 20.3 

0J4. 13.0 21.4 
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